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\ 7‘The  second  routine  augments  the  first  by  spelling  out  the  combined  effect  of  solar  heating  and 
solar  mass  attractiomjThe  semidiurnal  term  as  well  as  the  diurnal  term  are  considered. 

^The  third  algorithm  accounts  for  the  ocean  tide.  Only  the  semidiurnal  lunar  tide  is  presently 
considered,  but  the  algorithm  sets  the  pattern  for  introducing  various  other  spectral  components  of 
the  fluid  tide,  if  desired.  The  algorithm  starts  from  a table  for  the  global  tide  amplitudes  and  phase 
angles.  It  subsequently  approximates  the  gravitational  potential  exterior  to  the  earth  caused  by  the 
tidal  bulge  by  that  of  a system  of  point  masses  located  on  the  ocean  surface. 

The  fourth  tide  routine  is  optional;  Xx  may  be  used  in  place  of  the  third  algorithm  if  one  wishes 
to  base  the  tide  potential  on  a functionalization  of  the  tidal  sea  surface. 


Finally,  there  is  a computer  routine  for  the  Newtonian  attraction  caused  by  the  tide  bulge  of  the 
solid  earth.  Both  the  lunar  and  solar  tide  components  are  included.  The  tidal  properties  are  assumed 
to  be  functions  of  latitude.  Accordingly,  Love  coefficients  are  featured  that  are  zonal  harmonic 
expansions  of  latitude.  -Arf^llowance  has  been  made  for  tidal  lag.  The  latter  manifests  itself  as  a 
time  delay  in  the  response  of  the  tidal  mass  redistribution  to  the  tidal  stress  field.  Its  origin  are  the 
elasticity  and  plasticity  o\  the  earth  and  the  inertia  of  the  masses  involved. 
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FOREWORD 


Recent  progress  in  instrumentation  for  data  acquisition  and  computing  methods  suggests 
that  satellite  geodesy  must  concern  itself  with  even  such  minute  effects  as  the  gravitational  j 

action  that  the  tidal  dislocations  of  the  earth’s  masses  exert  upon  an  artificial  satellite.  A 
number  of  computer  algorithms  suited  to  this  purpose  were  developed  at  the  Naval  Surface 
Weapons  Center  (NSWC)  during  the  past  few  years.  The  relevant  physics  background  as  well 
as  sketches  of  the  computer  routines  were  published  in  a sequence  of  four  technical  reports 
in  1976  and  1977.  Since  then,  the  algorithms  were  worked  out  in  detail,  coded  for  the  • 

computer,  and  subjected  to  the  customary  numerical  checkout.  Also,  an  alternate  procedure  > 

was  developed  for  the  ocean  tides,  which  increased  the  number  of  tide,  routines  from  the 

original  four  to  five.  The  necessary  derivations  for  the  additional  routine  were  never 
documented.  They  are  included  as  Appendix  A.  The  present  report  is  intended  to  serve  as 
a manual  of  the  five  algorithms,  each  being  presented  in  a form  ready  for  immediate 
computer  implementation,  complete  with  checkout  data. 

The  work  documented  here  was  done  in  the  Astronautics  and  Geodesy  Division  and 
was  funded  as  part  of  the  development  of  the  TERRA  computer  program  which  is  used  to 

determine  the  earth’s  gravity  field  from  terrestrial  and  satellite  observations.  Henry  E. 

Castro  of  the  Computer  Program  Division  coded  those  algorithms  that  involve  the 
atmospheric  tides  and  the  two  versions  of  the  ocean  tide.  R.  Gordon  Barker,  also  of  the 
Computer  Program  Division,  coded  the  solid  earth  tide.  Additionally,  these  two  gentlemen 
invested  a great  effort  while  making  the  required  checkout  runs  on  the  computer.  Raymond 
B.  Manrique  of  the  Astronautics  and  Geodesy  Division  and  the  author  jointly  outlined  the 
computer  algorithms,  with  the  author  bearing  the  overall  responsibility  for  the  finished 
formulations  and  the  sole  responsibility  for  the  numerical  checkout  calculations. 

This  report  was  reviewed  by  Richard  J.  Anderle,  Head,  Astronautics  and  Geodesy 
Division. 


Released  by: 


R.  A.  NIEMANN,  Head 
Strategic  Systems  Department 


CONTENTS 


Page 

INTRODUCTION  1 

1 -COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION  CAUSED 

BY  THE  LUNAR  AIR  TIDE  2 

INPUT  DATA  2 

ALGORITHM  3 

COMPUTER  PROGRAM  OUTPUT 6 

TRIAL  DATA  FOR  PROGRAM  CHECKOUT 6 

2 - COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION  CAUSED 

BY  THE  SOLAR  AIR  TIDE 9 

INPUT  DATA  9 

ALGORITHM  10 

COMPUTER  PROGRAM  OUTPUT 14 

TRIAL  DATA  FOR  PROGRAM  CHECKOUT 14 

3 - FIRST  COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION 

CAUSED  BY  THE  M2  COMPONENT  OF  THE  OCEAN  TIDE  (BASED 

ON  A DISCRETE  REPRESENTATION  OF  THE  SEA  SURFACE) \7 

INPUT  DATA  17 

SUBROUTINE  FOR  THE  TIDE  POTENTIAL  COEFFICIENTS  18 

SUBROUTINE  FOR  THE  PERTURBING  ACCELERATION 23 

COMPUTER  PROGRAM  OUTPUT 28 

TRIAL  DATA  FOR  PROGRAM  CHECKOUT 28 

4 - SECOND  COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION 

CAUSED  BY  THE  M2  COMPONENT  OF  THE  OCEAN  TIDE  (BASED 

ON  A FUNCTIONALIZATION  OF  THE  SEA  SURFACE)  31 

INPUT  DATA  31 

SUBROUTINE  FOR  THE  TIDE  POTENTIAL  COEFFICIENTS  32 

SUBROUTINE  FOR  THE  PERTURBING  ACCELERATION 33 

COMPUTER  PROGRAM  OUTPUT 33 

TRIAL  DATA  FOR  PROGRAM  CHECKOUT 33 

5 - COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION 

CAUSED  BY  THE  TIDES  OF  THE  SOLID  EARTH  39 

INPUT  DATA  39 

SUBROUTINE  FOR  THE  POSITION  OF  THE  FICTITIOUS  MOON 41 

SUBROUTINE  FOR  THE  POSITION  OF  THE  FICTITIOUS  SUN  43 

SUBROUTINE  FOR  SATELLITE  POSITION  43 

V 


A 


CONTENTS  (Continued) 


Page 

SUBROUTINE  FOR  THE  PERTURBING  ACCELERATION 44 

COMPUTER  PROGRAM  OUTPUT 50 

TRIAL  DATA  FOR  PROGRAM  CHECKOUT 50 

DEVELOPMENT  STATUS  OF  THE  COMPUTER  ROUTINES  57 

REFERENCES 58 

APPENDICES 

A -NOTES  CONCERNING  THE  SECOND  OCEAN  TIDE  ROUTINE  61 

B - A COMPARISON  OF  THE  PERTURBING  ACCELERATION  DUE  TO 
THE  SOLID  EARTH  TIDE  WITH  THAT  CAUSED  BY  THE 
OCEAN  TIDE 69 

DISTRIBUTION 


vi 


INTRODUCTION 


Several  years  ago,  while  compiling  a force  model  on  which  to  base  the  trajectory 
integration  in  the  TERRA1  and  CELEST2  computer  programs  for  satellite  geodesy  and 
satellite  ephemeris  computation,  the  author  faced  the  nearly  complete  lack  of  published 
algorithms  that  would  in  a suitable  form  express  the  gravitational  action  by  which  the  tidal 
dislocation  of  the  earth’s  masses  perturbs  satellite  motion.  At  the  same  time,  it  was  realized 
that  the  earth  tides  would  have  to  be  taken  into  account  if  TERRA  and  CELEST  were  to 
substantially  exceed  the  accuracy  requirements  of  their  forerunners,  GEO3  and  ASTRO.4 
An  effort  was  consequently  made  to  derive  the  necessary  disturbing  potentials  and 
perturbing  accelerations.  The  result  consisted  of  a number  of  mathematical  models,  each 
designed  to  be  realistic  and  describing  one  type  of  perturbation  corresponding  to  one  of 
the  four  earth  tides.  The  details  of  these  models  were  documented  in  four  technical 
papers.5 ,6,7,8  Simultaneously,  the  necessary  computer  routines  were  developed  and 

submitted  for  coding. 

The  first  tide  routine  formulated  the  perturbing  acceleration  which  the  semidiurnal 

lunar  component  of  the  atmospheric  tide  bulge  exerts  on  the  orbit  of  an  artifical  satellite. 
The  tide  bulge  was  assumed  to  result  from  the  fact  that  the  earth  rotates  within  the  field 

of  lunar  mass  attraction  which  is  inhomogeneous  across  the  terrestrial  globe.  The  routine 

was  restricted  to  the  main  term  of  the  semidiurnal  tide. 

The  second  tide  routine  was  designed  to  augment  the  first  by  spelling  out  the 
combined  effect  on  the  atmosphere  of  solar  heating  and  solar  mass  attraction,  plus  the 
resulting  modification  of  the  earth’s  gravitation.  The  semidiurnal  term  as  well  as  the  diurnal 
term  were  included. 

The  third  algorithm  was  to  account  for  the  ocean  tide.  It  started  from  an  existing 
table  for  the  global  tide  amplitudes  and  phase  angles,  and  postulated  that  each  tidal  height 
value  resembles  a gravitating  point  mass  that  by  its  presence  perturbs  the  satellite  motion. 
Only  the  semidiurnal  lunar  tide  was  considered;  however,  this  algorithm  may  be  regarded  as 
setting  the  pattern  for  introducing  various  other  spectral  components  of  the  fluid  tide,  if 
desired. 

Further,  there  was  a computer  routine  for  the  Newtonian  attraction  caused  by  the  tide 
bulge  of  the  solid  earth.  Both  the  lunar  and  the  solar  tide  components  were  included.  The 
tidal  properties  were  assumed  to  be  functions  of  latitude.  Accordingly,  Love  coefficients 
were  featured  that  are  zonal  harmonic  expansions  of  latitude.  An  allowance  was  made  for 
tidal  lag  manifesting  itself  as  a time  delay  in  the  response  of  the  tidal  mass  redistribution 
to  the  tidal  stress  field. 
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Several  months  after  completing  the  coding  requests  for  these  four  algorithms,  an 
opportunity  presented  itself  to  devise  an  interesting  alternative  to  the  just  mentioned  ocean 
tide  routine.  Namely,  a functionalization  became  available  for  the  tidal  ocean  surface.  This 
consisted  of  a representation  of  the  tidal  ocean  surface  by  an  expansion,  in  surface 
harmonics,  that  had  been  least  squares  fit  adjusted  to  the  table  of  discrete  tide  amplitudes 
and  phases  that  had  been  the  basis  of  the  ocean  tide  algorithm.  Upon  receiving  a card 
deck  containing  the  coefficients  for  this  expansion  [made  available  to  us  through  the 
kindness  of  Dr.  C.  C.  Goad  of  the  National  Oceanic  and  Atmospheric  Administration 
(NOAA)] , an  optional  routine  was  written,  which  interprets  the  functionalized  ocean  surface 
as  a continuous,  gravitating,  tidal  mass  layer.  From  it,  the  perturbing  acceleration  acting  on 
the  satellite,  due  to  the  presence  of  the  oceanic  tide  bulge,  is  then  derived  via  the  Poisson 
integral  and  subsequent  gradient  formation.  While  the  physics  background  and  the 
mathematical  derivations  for  the  original  four  tide  routines  are  available  in  the  shape  of 
formal  reports,5 ,6, 7,8  the  background  for  this  latter  optional  version  of  the  ocean  tide 
algorithm  had  never  before  been  documented.  Thus,  to  preserve  it,  it  is  included  in  this 
report  as  Appendix  A. 


1 - COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION 
CAUSED  BY  THE  LUNAR  AIR  TIDE 


INPUT  DATA 

R = any  reasonable  value  for  the  mean  radius  of  the  earth,  in  meters 
G = gravitational  constant,  in  meters,  kilograms,  and  seconds 
suggested  value:  G = 6.6732E  - 11  m3  kg'*  sec-2 
A2  = a constant  associated  with  the  physics  of  the  lunar  air  tide 
suggested  value:  A2  = 0.564  kg  nT2 

The  parameters  J,  n,  and  t*  specify  the  time  instant  for  which  the  perturbing  acceleration 
due  to  the  tide  is  to  be  found.  Normally,  this  is  the  time  instant  associated  with  the 
current  time  line  of  the  orbit  integration. 

J = number  of  the  calendar  year 

n = number  of  the  day  within  the  year 

t*  = mean  solar  time  at  Greenwich  (GMT,  UTC),  in  seconds 

x,|  inertial,  Cartesian  components  of  satellite  position  vector,  in  meters, 
x2>=  associated  with  the  just  specified  time  instant  (for  a precise  definition  of 
x3J  the  reference  frame,  see  Appendices  B and  C of  Reference  2) 
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ALGORITHM 


N 


6(ij)  = 


if 


= J 


(1-1) 


10  i 

= n 6(J,  1975)  + (365  + n)8(J,  1976)  + (731  + n)  8(J,  1977)  + (1096  + n)  6(J,  1978) 


+ (1461 +n)5(J,  1979)  + (1826 + n)8(J,  1980)  (1-2) 

AT  = 5.28E  - 04  + (3.56E  - 08)N  (1-3) 

d = 27392.5  +N  + AT  + t*  (M) 

T = d/36525  (1-5) 

s = 270.434358  + 48 1267.883 14  IT-  0.001 133T2  + 0.000002T3-  (1-6) 


h = 279.69668  + 36000.768930T  + 0.000303T2 


s and  h will  result  in  decimal  fractions  of  degrees. 

Calculate  now  the  earth-fixed,  Cartesian  satellite  coordinates  (y15  y2,  and  y3)  from  the 
corresponding  inertial  coordinates  (x,,  x2,  and  x? ). 


= (ABCD)j  n ,t* 


(1-8) 


where  (ABCD)j  n is  the  transformation  matrix  that  manages  the  transition  from  the 

inertial  frame  of  the  satellite  equations  of  motion  ( Basic  Inertial  System,  which  is 

associated  with  either  1950.0  or  with  the  beginning  of  the  day  UT  of  trajectory  epoch)  to 
the  earth-fixed  reference  frame  corresponding  to  the  time  instant,  J,  n,  t*.  This 
transformation  is  a computer  routine  that  is  already  part  of  TERRA/CELEST.  It  is  external 
to  the  present  algorithm.  For  details,  see  Appendices  B and  C of  Reference  2. 
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In  the  earth-fixed  frame,  the  components  of  the  perturbing  acceleration  are 


where 


3U 

3r 

3U 

3X 

3U 

30 

3r 

9yj 

' 3X 

3yi 

' 30 

9yi 

a> 

C 

3r 

3U 

3X 

3U 

30 

dr 

3y2 

' 3X 

3y2 

' 30 

9y2 

3U 

3r 

3U 

3X 

3U 

30 

3r 

dy3 

' 3X 

9y3 

' 30 

9y3 

(1-9) 


(1-10) 


(Ml) 


3r 

= - i = 

r 

3X 

y2 

9yi 

y?  + y\ 

3X 

yi 

dy2 

' y\  +y2 

3X 

9y3 

= 0 

30 

-yi  yj 

3y, 

r2  vV2  + y2 

30 

- y2  yj 

dy2 

r2  y/y]  + y\ 

30 

dy3 

s/y\  + y2 
r2 

(1-12) 

(1-13) 

(1-14) 

(1-15) 

(1-16) 

(1-17) 

(1-18) 
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(1-19) 


9U  3R3  , 5R5  , 

= - a ~p-  Pj  (sin  0)  cos  2a  + b -p-  P2(sin  0)  cos  2a 

dU  2R3  , 2R5  , 

^ = - a -jy  P^sin  0)  sin  2a  + b -p-  P2(sin  0)  sin  2a 

dU  3r3  . 1 5 R5  , 

^ ~ ~ a —f  sin  20  cos  2a  + b -p—  (7  sin2  0-4)  sin  20  cos  2a 


r = 


a = 


b = 


+ ^ y]  + y\  + y\ 

Sir2 

A,GR 

2 64 

5n2 

A,GR 


3072  48 


P2,(x) 


3(1  - x2) 
15 


P4M  = — (1  - x2)(7x2  - 1) 


x = 


• n y3 
sin  0 = — 


sin  20 


2y3Vv  1 +V2 


cos  2a 
sin  2a 


cos2  a - sin2  a 


2 sin  a cos  a 
cos  a = cos  a*  cos  X - sina*  sin  X 
sin  a = sin  a*  cos  X + cos  a*  sin  X 


a’  = 


v = 


cos  X = 


t**-  y-7?5 
360 

86400  1 
s - h 

y2 


Vyf  + 


(1-20) 

(1-21) 

(1-22) 

(1-23) 

(1-24) 

(1-25) 

(1-26) 

(1-27) 

(1-28) 

(1-29) 

(1-30) 

(1-31) 

(1-32) 

(1-33) 

(1-34) 

(1-35) 

(1-36) 
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sin  X 


_Zi_ 
Vy2  + y\ 


Finally,  transform  the  perturbing  acceleration  back  into  the  inertial 


where  (ABCD)T  is  the  transpose  of  (ABCD). 


COMPUTER  PROGRAM  OUTPUT 

The  inertial  components  (Tx  ] , TX2>  and  TXj)  of  the  perturbing 
the  presence  of  the  lunar  air  tide  are  output  in  terms  of  m sec"2.  If 
km  sec"2. 


TRIAL  DATA  FOR  PROGRAM  CHECKOUT 

R = 6378145  m 

G = 6.6732E  - 1 1 m3  kg" 1 

A2  = 0.564  kg  m-2 

J = 77 

n = 202 

t*  = 50000  sec 

x,  = 3151529.23  m 

x2  = 5458608.75  m 

x3  = 3639072.50  m 


(1-37) 


frame 


(1-38) 


acceleration  due  to 
required,  convert  to 


sec 
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From  external  source  (CELEST  subroutine  PRENUT) 

/-  0.8405285753  0.5417623775  0.2289080162E  - 02 

(ABCD)  = - 0.5417605355  -0.8405316908  0. 1413662999E  - 02 

\ 0.26899 13850E-  02  -0.5 190827376E-  04  0.9999963803 

y,  = +316648.61  m 

y2  = - 6290363.38  m 

y3  = +3647253.32  m 

t**  = 208?3333333 

N = (73 1 + 202)  = 933 

AT  = 5.612148000E-  04 

t*  = 0.5787037037  deg 

d = 2.832607926E  + 04 
T = 0.7755257840 
h = 2?8 19922 141 E + 04 
s = 3?73506086 1 E + 05 
X = 272?88 1 76 1 6 
h - s = - 3?4530686469E  + 05 
a = -3?448331496E  + 05 
2a  = - 6?896662992E  + 05 
x = 0.5011240256 
P^(x)  = 2.246624133 
P^(x)  = 4.256662025 
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au 

3r 

= 7.0700492 16E-  12 

m sec’ 

au 

ax 

= - 5.416307419E  - 04 

m2  sec 

au 

30 

= 2.46749468 IE  - 05 

m2  sec 

3r 

9y, 

= 4.350677408E  - 02 

ax 

dyj 

= 1.5857151 04 E - 07 

m’1 

30 

ay, 

= - 3.4615995  18E-  09 

nT1 

3r 

3y2 

= - 8.64281 1299E-  01 

ax 

ay2 

= 7.98228 1040E-  09 

m'1 

30 

ay2 

= 6.8766191 14E  - 08 

nT1 

3r 

9y3 

= 5.01 1240258E  - 01 

ax 

9y3 

= 0 

30 

ay3 

= 1 . 189005543E  - 07 

nT  1 

Ty, 

= - 8.566502457E  - 11 

m sec-2 

TV2 

= - 8.737156821E-  12 

m sec'2 
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Ty  = 6.476836379E  - 12 


Tx  = 7.675476994E  - 1 1 
‘ = 7.675476994E  - 14 

Tx  2 = - 3.906656638E  - 11 
= - 3.906656638E  - 14 

TX3  = 6.26836743 1 E - 12 
= 6.268367431  E-  15 

T = |f|  = 8.635267078E  - 1 1 
= 8.635267078E  - 14 


m 

sec-2 

m 

- 2 
sec 

km 

- 2 
sec 

m 

sec"2 

km 

sec"" 

m 

sec  i 

km  sec" 

m 

- 2 
sec  * 

km  sec" 

2 - COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION 
CAUSED  BY  THE  SOLAR  AIR  TIDE 


INPUT  DATA 

R = any  reasonable  value  foT  the  mean  radius  of  the  earth,  in  meters 
G = gravitational  constant,  in  meters,  kilograms,  and  seconds 
r suggested  value:  G = 6.6732E  - 11  m3  kg-1  sec-2 

Aj  j constants  associated  with  the  physics  of  the  solar  air  tide 
) = suggested  values:  Aj  = 6 kg  m"2 
A2  j A2  = 11.9  kg  m'2 

The  parameters  J,  n,  and  t*  specify  the  time  instant  for  which  the  perturbing 
acceleration  due  to  the  tide  is  to  be  found.  Normally,  this  is  the  time  instant  associated 
with  the  current  time  line  of  the  orbit  integration. 

J = number  of  the  calendar  year 
n = number  of  the  day  within  the  year 
t*  = mean  solar  time  at  Greenwich  (GMT,  UTC),  in  seconds 

x,  i inertial,  Cartesian  components  of  the  satellite  position  vector,  in  meters, 
x2  7 = associated  with  the  just  specified  time  instant  (for  a precise  definition 
x3  \ of  the  reference  frame,  see  Appendices  B and  C of  Reference  2) 
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ALGORITHM 


r 


u 


First,  calculate  the  earth-fixed,  Cartesian  satellite  coordinates  (y,,  y2,  and  y3 ) from 
the  corresponding  inertial  coordinates  (xr  x2,  and  x3). 


= (ABCD>J,n,.*  (2-D 

where  (ABCD)j  n is  the  transformation  matrix  that  manages  the  transition  from  the 
inertial  frame  of  the  satellite  equations  of  motion  ( Basic  Inertial  System , which  is 
associated  with  either  1950.0  or  with  the  beginning  of  the  day  UT  of  trajectory  epoch)  to 
the  earth-fixed  reference  frame  corresponding  to  the  time  instant,  J,  n,  t*.  This 
transformation  is  a computer  routine  that  is  already  part  of  TERRA/CELEST.  It  is  external 
to  the  present  algorithm.  For  details,  see  Appendices  B anc  C of  Reference  2. 

In  the  earth-fixed  frame,  the  components  of  the  perturbing  acceleration  are 


T = ^ 3U  3X_  dU  30 

y‘  dr  dy,  3X  0y , 30 

au  dr + au  ax_  au  a o 

dr  dy2  3X  3y2  39  dy2 

X = 9U9r_  + 3UaiX_  3U  30_ 

J 3 3r  3y3  3X  3y3  36  3y3 


where 


3r_  = y, 
dyj  r 


i = 1,2,3 


ax_  = y2 

9yi  y]+y\ 

ax  = y, 
ay2  y?  +y2 


(2-2) 

(2-3) 

(2-4) 


(2-5) 

(2-6) 

(2-7) 
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ax 

0y3 


— = o 


dj_ 

ay, 

ae 

ay2 

be 

ay3 

au 

3r 


au 

ax 


au 

30 


- y.  y3 

rVyf  + y2 

- y3  y3 

rVy]  + y^ 

v'y  j + y2 
r2 

4R4  -JR3 

a,  -j3-Pj(sin0)  cos  a - a2  — j-  P2  (sin  0 ) cos  20 

5r5 

+ a3  — P2(sin  0)  cos  20 

R4  2R3 

ai  ^Pj(sin0)  sin  a - a2  — j-  P2(sin0)  sin  20 

2r5 

+ a3  — j-  P2(sin  0 ) sin  20 

* 

3R4  3r3 

a,  (15  sin2  0 - 11)  sin  0 cos  a - a2  -j—  sin  20  cos  20 


r = 


a,  = 


a-,  = 


a,  = 


Pi(x)  = 


(2-8) 


(2-9) 


(2-10) 


(2-11) 


(2-12) 


(2-13) 


+ a3  ^ (7  sin2  0-4)  sin  20  cos  20 

(2-14) 

= + v/y?  +y2  +y] 

(2-15) 

8ttGR 

= A, 

105 

(2-16) 

5tt2GR 

= 

2 64 

(2-17) 

5tt2GR  a2 

2 3072  " 48 

(2-18) 

= - vl  - x2(5x2  - 1) 

(2-19) 

1 1 


(2-20) 


P|(x)  = 3(1  - x2) 

P4 (x)  = y (1  - x2 )(7x2  - 1) 


(2-2 1 ) 


_ . a V3 
x = sin  a = — 


sin  28  = 


2y3Vy?  +y  \ 


(2-22) 


(2-23) 


cos  a = cos  a*  cos  X - sin  a*  sin  X 


(2-24) 


sin  a = cos  a*  sin  X + sin  a*  cos  X 


(2-25) 


cos  2 P = cos2  P - sin2  P 


(2-26) 


sin  2p  = 2 sin  P cos  P 


(2-27) 


cos  p = cos  P*  cos  X - sin  P*  sin  X 


(2-28) 


sin  P = sin  P*  cos  X + cos  P*  sin  X 


(2-29) 


cos  X = 


'y\  +y? 


(2-30) 


sin  X = 


y r^i 


(2-31) 


a*  = t**  - 78° 


(2-32) 


P * = t**-  146° 


(2-33) 


360 

t**  = 1* 

86400 


(2-34) 


Finally,  transform  the  perturbing  acceleration  back  into  the  inertial  frame 


1 \ t ( yi 

r 1 — 1 a nrniT  x 


TX2  j-  (ABCD)‘  n t.  Ty2 


(2-35) 
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where  (ABCD)T  is  the  transpose  of  (ABCD). 

While  coding  the  above  algorithm  for  use  with  TERRA  and  CELEST,  it  was  found 
helpful  to  keep  in  mind  the  following.  A time  line  is  characterized  by  the  time  elapsed 
since  the  epoch  (starting  instant)  of  the  particular  orbital  arc.  The  epoch  is  specified  in 
terms  of  year,  number  of  the  day  UT  in  the  year,  and  time  (tEP),  in  seconds,  elapsed 
from  the  start  of  the  day  UT  until  epoch  time.  Usually,  tEp  =0.  The  present  algorithm  is 
periodic  in  t*,  the  latter  being  the  Universal  time  associated  with  the  individual  time  line, 
t*  starts  from  zero  at  the  beginning  of  the  day  UT  during  which  the  time  line  occurs.  It 
is  thus  unimportant  how  many  integer  days  have  elapsed  since  epoch.  For  an  illustration, 
see  Figure  1. 


EPOCH 


*EP 


EPOCH 
DAY  UT 


At 


INTEGER  NUMBER  OF  DAYS 


TIME  LINE 


DAY  UT 
OF 

TIME  LINF 


Figure  1.  Universal  Time  in  the  Algorithm  for  the 
Perturbing  Acceleration  Due  to  the  Solar  Air  Tide 


To  actually  find  the  numerical  value  of  t*  corresponding  to  the  time  line,  consider 

that 


tEP  = Epoch  time  (Universal  time  at  epoch) 

At  = time  difference  between  time  line  and  epoch 
t*  = Universal  time  of  time  line 

Now,  express  tEp  (if  non-zero)  and  At  and  t*  in  terms  of  a common  time  unit  (usually 
mean  solar  seconds).  Find 
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t'  = tEp+At  (2-36) 

Express  t'  in  terms  of  days  and  decimal  fraction  of  days  (1  day  = 86400  sec).  Now, 

t*  = FRACT(t')  (2-37) 

Finally,  convert  t*  to  angular  time 

t**  = 360°  FRACT(t')  (2-38) 

which  may  be  substituted  for  Equation  2-34. 

COMPUTER  PROGRAM  OUTPUT 

The  inertial  components  (TXj,  TX2,  and  TXj)  °f  the  perturbing  acceleration  due  to 
the  presence  of  the  solar  air  tide  are  output  in  terms  of  m sec-2.  If  required,  convert  to 
km  sec-  2 . 

TRIAL  DATA  FOR  PROGRAM  CHECKOUT 

R = 6378145  m 

G = 6.6732E  - 11  m3  kg-1  sec-2 

A,  = 6 kg  m-2 

A2  = 1 1.9  kg  m-2 

J = 1977 

n = 202 

t*  = 50000  sec 


14 


X1 

= 3151529.23 

m 

X2 

= 5458608.75 

m 

= 3639072.50 

m 

From  external  source  (CELEST  subroutine  PRENUT) 


(ABCD)  = 


0.8405285753 
0.5417605355 
0.26899 13850E  - 02 


0.5417623775 
“0.8405316908 
“0.5190827376E-  04 


0.2289080 162E  - 02 
0.1413662999E  - 02 
0.9999963803 


Vl 

y2 


y3 


X 

9 


t** 


0 

P|(x) 

P\(x) 

P^(x) 


316648.61  m 

- 6290363.38  m 

3647253.32  m 

7278145.00  m 

272°8817616 

30°07439285 

208*13333333 

403° 2 15095 

335°2 15095 


0.3318192840 

2.246624133 

4.256662025 


6.1 12661413E  - 04 

m2 

sec-2 

3.905397704E  - 03 

m2 

sec' 2 

8. 1 362452 1 7E  - 05 

m2 

sec' 2 

15 


— = - 1.450823033E-  09 
dr 


>2 


‘y-< 


2 


T 


X 


3 


T 


8.662262286E  - 10 

m sec-2 

- 1 .5 188275 19E  - 09 

m sec-2 

- 1.612467289E-  09 

m sec-2 

- 1 .6 12467289E  - 12 

km  sec- 

6.1912321 15E-  12 

m sec-2 

6.1912321 15E  - 15 

km  sec- 

- 1 .5 14495280E  - 09 

m sec-2 

- 1 .5 14495280E  - 12 

km  sec- 

1T|  = 2.2 1 2 1 90 1 0 1 E - 09 

m sec- 1 

= 2.2121 90 1 0 1 E - 12 

km  sec- 

3 - FIRST  COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION 
CAUSED  BY  THE  M2  COMPONENT  OF  THE  OCEAN  TIDE  (BASED  ON  A 
DISCRETE  REPRESENTATION  OF  THE  SEA  SURFACE) 


INPUT  DATA 

= tidal  amplitude  on  area  element  ij  (these  are  output  data  from  the 
NSWC/Schwiderski  Ocean  Tide  Program);  if  available  in  meters,  use 
directly;  otherwise,  convert  to  meters 

= tidal  phase  angle  on  area  element  ij  (these  are  output  data  from  the 
NSWC/Schwiderski  Ocean  Tide  Program),  available  in  terms  of  degrees 

NP  = number  of  grid  points  in  the  NSWC/Schwiderski  Ocean  Tide  Program; 
this  may  be  expected  to  be  a five-digit  integer 

R = radius  of  the  earth  (semimajor  axis  of  a suitable  reference  ellipsoid), 
in  kilometers 

suggested  value:  R = 6378.145  km 
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e2  = square  of  eccentricity  of  reference  ellipsoid 
suggested  value:  e2  = 0.00669342 

in  case  it  is  desired  to  start  from  the  flattening,  f,  find  e2  from 
e2  = (2  - Of 

Me  = gravitational  constant  of  the  earth,  in  kilometers  and  seconds 
suggested  value:  jufc  = 398601  km3  sec-2 

G = Newton’s  constant,  in  kilometers,  kilograms,  and  seconds 
suggested  value:  G = 6.6732E-20  km3  kg-1  sec-2 

p = density  of  water,  in  kilograms  and  kilometers 
suggested  value:  p = l.E+12  kg  km-3 

a - rate  of  mean  mean  longitude  of  moon 

180  . 
suggested  value:  a = — 1.405 19E-  04  deg  sec  ' 

7T 

The  parameters  J,  n,  and  t*  specify  the  time  instant  for  which  the  perturbing 
acceleration  due  to  the  tide  is  to  be  found.  Normally,  this  is  the  time  instant  associated 
with  the  current  time  line  of  the  orbit  integration 

J = number  of  the  calendar  year 
n = number  of  the  day  within  the  year 
t*  = mean  solar  time  at  Greenwich  (GMT,  UTC),  in  seconds 


x,\  inertial,  Cartesian  components  of  satellite  position  vector,  in 
x2  ( = kilometers,  associated  with  the  just  specified  time  instant  (for  a 
x3  ( precise  definition  of  the  reference  frame,  see  Appendices  B and  C of 
' Reference  2) 


NMAX  = limit  on  the  number  of  terms  in  the  expansion  of  the  tide  potential 


SUBROUTINE  FOR  THE  TIDE  POTENTIAL  COEFFICIENTS 

This  is  the  first  part  of  the  algorithm  for  the  first  ocean  tide  computer  routine.  The 
time-independent  constituents  of  the  expansion  coefficients  of  the  ocean  tide  potential 
(“Fnm>  0Fnm’  ^nm’  and  are  calculated  from  the  amplitude  and  phase  tables 

produced  by  the  NSWC/Schwiderski  Ocean  Tide  Computer  Program.  Once  established,  the 
latter  remain  unchanged  as  long  as  the  current  versions  of  the  amplitude  and  phase  tables 
remain  valid.  Thus,  the  subroutine  under  discussion  will  be  exercised  quite  infrequently, 
namely,  just  once  each  time  a new  edition  of  the  just  mentioned  amplitude  and  phase 
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tables  becomes  available.  That  is  expected  to  occur  at  intervals  of  several  years,  during 
which  the  present  subroutine  will  remain  dormant. 

For  the  geometry  associated  with  the  index  ij,  see  Figures  2,  3,  and  4. 


a{j  = 10~3p  G ASjj 

fit 

cos  6jj 

(3-1) 

0ij  = 10-3p  G AS^ 

fit 

sin  8y 

(3-2) 

/ vr  \2  2 

/ 7T  \ 

'>  ~ (iso)  R 

sin 

oc  I 
0 1 

(3-3) 

= area  of  the  nonpolar  surface  area  element  for 
j = 2,3,4  •••  ,jMAX<180 

3 


(ASuWar  2 ( 1 80 ) R 


(3-4) 


= area  of  the  polar  surface  area  element  (neighboring  the  North 
Pole) 


„ n n l 1 \ 

N = TiSV'  2)  “,r“'lianS 

N ■ Fso  (‘  ■ 2)  ln  rad,ans 

sp"J*u) 


(3-5) 

(3-6) 

(3-7) 


To  each  index  ij  assign  now  a counting  number  v,  1 < v < NP. 


'H 


(n  - in)!  1 


1 

My  R2  " 

N P 

T p2n+la  P 

v = 1 

(3-8) 

1 

pER2n 

N P „ 

r p2n+l0  p 

v-  1 

(3-9) 

NP 

y p2n+1a  h" 

v nm 

v ~ 1 

(3-10) 
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y,  INDICATES  THE  POINT  WHERE  GREENWICH  MERIDIAN  AND  EQUATOR  INTERSECT 
y2  ' LOCATED  IN  THE  EQUATORIAL  PLANE 


Figure  4.  Position  of  the  Point  Mass  M,j  in  the 
Earth-Fixed  Cartesian  Coordinate  Frame 


_ Mil 

(n  + m)!  MER2n  fe\  v u nr 


K = 


if 


B = k 
t * k 


(3-11) 


(3-12) 


Calculate  the  f*V  and  h*  as  follows,  noting  that 

n m n m ° 


R 


P»  = 


(3-13) 


and 


hv  =0  for  all  values  of  n 

n o 


(3-14) 
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r 


and  that 


► 

t 


L 


li 


V 

n m 


= r = o 

n m 


for  any  n < m 


(3-15) 


Obtain,  separately  for  each  v , the  required  f"m  and  h*m  from  the  following 
recurrence  relations: 

To  advance  in  n,  evaluate 


f"*  ■ .-  * [<2n  * 1 ■ > sin  °.  - 1-  * ->  p.  c. , . J 6> 

h»*l,m  = n-  m+  I [<2"  + U Sl"  "•  ' ("  + m)  P-  h»-  ' J (3'l7) 

To  advance  in  m,  use 

fn+ 1 ,n+ 1 = (2n  + l)p„  [cos  6v  cos\v  F n - cos  6v  sin  X^  h’J  (3-18) 

hn+ 1 ,n+ 1 = (2n+Dpi([cos0vcosXi;h;jn+cos0vsinXi,Pjn]  (3-19) 

Start  these  recurrences  from 


‘o.o 


p 

1,0 


= R 


sin  0„ 


(3-20) 

(3-21) 


h 


V 

0,0 


hP0  = 0 


(3-22) 


and  terminate  the  procedure  at  n = NMAX. 


Store  the  resulting  aFnm,  (,Fnm,  aHnm,  and  0Hnm . They  will  be  programmed  as 
constant  parameters  into  the  subroutine  for  the  perturbing  acceleration  and  will  be  expected 
to  serve  for  all  subsequent  runs  of  that  subroutine,  until  updated.  Note  that  the  “F 
^nm'  a Hn m ’ anc*  ^nm  are  dimensionless  quantities. 
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While  it  would  be  decidedly  impractical  to  compute  f"m  and  h"m  from  the  following 
equation,  for  program  checkout  purposes  it  is  possible  to  find  individual,  numerical  values 
for  these  quantities  from 


P™  (sin  0v) 


This  relationship  may  also  be  used  to  justify  Equation  3-15,  considering  that 

P™(x)  = 0 for  n < m 

which  is  one  of  the  basic  identities  of  Legendre  function  theory. 


(3-23) 


(3-24) 


SUBROUTINE  FOR  THE  PERTURBING  ACCELERATION 

First,  calculate  the  earth-fixed,  Cartesian  satellite  coordinates  (y , , y2,  and  y3 ) from 
the  corresponding  inertial  coordinates  (x,,  x2,  and  x3). 


= (ABCD)j  n t. 


(3-51)+ 


where  (ABCD)j  n t,  is  the  transformation  matrix  that  manages  the  transition  from  the 
inertial  frame  of  the  satellite  equations  of  motion  (Basic  Inertial  System,  which  is 
associated  with  either  1950.0  or  with  the  beginning  of  the  day  UT  of  trajectory  epoch)  to 
the  earth-fixed  reference  frame  corresponding  to  the  time  instant,  J,  n,  t*.  This 
transformation  is  a computer  routine  that  is  already  part  of  TERRA/CELEST.  It  is  external 
to  the  present  algorithm.  For  details,  see  Appendices  B and  C of  Reference  2. 

Find 


r = + (3-52a) 

or 

r = +Vypryf~+"yJ  (3-52b) 


as  convenient 


t 


Numbering  sequence  advanced  to  emphasize  parallelism  between  certain 


equations  of  the  algorithms  in  this  and  in  the  next  chapter. 
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Calculate  now  the  time  dependent  expansion  coefficients  for  the  tide  potential  from 


Fnm  = “Fnm  COS<at*  +X)  + (iFnm  Sinlot*  +X)  (3-53) 

Hnm  = a Hnm  COS  K + X>+<,Hnm  sin  (at*  + X)  (3*54) 

Obtain  x as  follows 

{1  n = m 

if  (3-55) 

0 n ¥=  m 

Nv  = N 5(J,  1975)  + (365  + N)  6(J,  1976)  + (73 1 + N)  S(J,  1977) 

+ (1096  + N)  6(J,  1978)  + (1461  + N)  6(J,  1979)  (3-56) 

AT  = 5.28E-  04  + (3.56E- 08)N£ 

d0  = 27392.5  + Nl  + AT 

X - do 
0 36525 

X = 270.434358  + 48  1267. 883 141 37  T„- .00 1133  T* +.00000 19  T;*  (3-60) 

Note  that  the  day  number  N,  and  thus  x,  must  be  updated  whenever  t*  runs  into  the 
next  day. 

Remember  that  Fnm  and  Hnm  are  linear  functions  of  two  trigonometric  terms.  To 
evaluate  these  trigonometric  terms,  let  t*  and  t*+ , be  the  values  of  Universal  time  for 
which  subsequent  integration  steps  (time  lines)  are  to  be  performed.  Update  the 
trigonometric  time  factors  as  follows: 

cos  (<Jt*+ , + x)  = cos  { (ot*  + x)  + oAt  } = cos  (at*  + x)  cos  aAt  - sin  (at*  + x)  sin  oAt  (3-6 1 ) 


(3-57) 

(3-58) 

(3-59) 
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sin  (ot*+ , + x)  = sin  { (at*  + x)  + oAt  } = sin  (at*  + x)  cos  oAt  + cos  (at*  + x)  sin  aAt  (3-62) 


**  - c.-c 


(3-63) 


Now  calculate  the  Cartesian  components  of  the  perturbing  acceleration  in  the 
earth-fixed  frame 


d6  N M A X n / 3U  3Y 

dy.  n=0  V nm  nm  3y, 


(3-64) 


N M A X n / 3U  3Y 

y y F — ^ + H — ^ 

i—i  V nm  av  nm  a 

n = 0 m =0  \ °’2  "V? 


(3-65) 


where 


NM  a x n / 

E E (Fn 

n = 0 m =0  \ 


m ay. 


nm  3y, 


1/1  1_ 

R \2  ^h,n^n+l.m-t  2 ^n  + 1 .m  + 1 


1 / ! 1 
— I _ - A V -—V 

R \ 2 mn  n+1>m-l  2 n+ 1 ,m  + 1 


— (n- m + l)Un+1  m 

1 (\  1 

- I - A V - - V 

£ ^2  m n n+l,m-l  n+l,m+l 


R \ 2 n ^n+  l»m-l  2 ^n+  1 .m  + 1 


— (n-  m + l)Vn+1  r 


(3-66) 


(3-67) 


(3-68) 


(3-69) 


(3-70) 


(3-71) 


(3-72) 
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1 

] 

I 

and 


U 


n 

= (n  - m + l)(n  - 

m + 2) 

(3-73) 

Un,-. 

(n  - 1 )! 

= - — — U , 
(n  + 1 )!  nl 

U„i 

n(n+  1) 

(3-74) 

Vn,-  1 

(n-  1)! 

= - — V 

(n  + 1 )!  nl 

V„1 

n(n  + 1 ) 

(3-75) 

Obtain  the  values  of  the  individual  eigenfunctions  from  the  following  recursive 
relationships: 


P 

R 

r 

(3-76) 

sin 

r 

(3-77) 

Vno 

= 0 for  all 

values  of  n 

(3-78) 

unm 

= Vnm  =0 

n m 

for  all  n < m 

(3-79) 

To  advance  in  n,  evaluate 

Un+l.n,  = n - m + 1 j(2n+  1)S'n^Unm  - (n  + m)PUn-.,ml  (3'80) 

Vn  + 1 ,m  = n_“  + 1 |(2n+  I)sin*  V nin  - (n  + m)p  Vn_  , m [ (3-81) 

To  advance  in  m,  use 

' «"+'>p(7  7V»n)  W'82) 


26 


V = 

v n+  1 ,n+  1 

(2n+  l)p(~  V + - Un) 

\ r nn  r nn  j 

(3-83) 

Start  from 

c 

O 

© 

II 

Me 

r 

(3-84) 

U.,o  = 

R y3 

Me  r3 

(3-85) 

II 

<o 

o 

> 

o 

n 

o 

(3-86) 

The  Cartesian  components  of  the  perturbing  acceleration,  in 

the  earth-fixed  reference 

frame,  are  now 

Tv,  = 

30 

(3-87) 

y i 

3yj 

T 

1 y 2 

00 

dy2 

(3-88) 

Tv  = 

00 

(3-89) 

y 3 

9y3 

Finally,  rotate  the  perturbing  acceleration  back  into  the  inertial  coordinate  system 


(3-90) 


where  (ABCD)T  is  the  transpose  of  (ABCD). 
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COMPUTER  PROGRAM  OUTPUT 

The  inertial  components  (Tx  ( , TXl,  and  Tx , ) of  the  perturbing  acceleration  due  to 
the  presence  of  the  ocean  tide  are  output  in  terms  of  km  sec  2.  If  required,  convert  to 
m sec*2. 


TRIAL  DATA  FOR  PROGRAM  CHECKOUT 

Non-zero  tidal  heights  and  associated  phase  values  were  assumed  for  nine  surface  area  • 
elements.  Several  of  the  “F  , SF  , “H„„  , and  were  subsequently  computed.  The 

algorithms  for  the  perturbing  acceleration,  featured  in  the  first  and  second  computer 
routines  for  the  ocean  tide,  are  identical.  Trial  data  for  these  segments  of  the  two  tide 
routines  will  therefore  be  listed  only  once,  namely,  for  the  second  ocean  tide  routine, 
below. 


Quantities  Associated  With  the  Surface  Area  Elements 


f 

6 

0 

X 

p 

AS 

V 

a 

V 

K 

ki 

V 

m 

deg 

rad 

rad 

km 

km2 

km1  sec-’ 

kmJ  sec'1 

t.i 

1 

10 

25 

1.562069681 

8.7266463  t 03 

6356.800824 

108.1411251 

6.5403462E  - 08 

3.0498135E  - 08 

2.1 

2 

10 

25 

1.562069681 

2.617993881.  02 

6356.800824 

108.141  1251 

6.5403462E  - 08 

3 0498135E  - 08 

3.1 

3 

10 

25 

1.562069681 

4.363323 13I£  - 02 

6356.800824 

108. 1411251 

6 5403462E  - 08 

3.0498 1 35E  - 08 

1.2 

4 

10 

25 

1.544616388 

8.7266463  E - 03 

6356  813825 

432  4766612 

2.6 1 56072E  - 07 

1.21 96  7 7 7 E - 07 

2.2 

5 

20 

30 

1 544616388 

2.6 1 799388E  - 02 

6356.813825 

432.4766612 

4.9987043E  - 07 

2.8860033E  - 07 

3.2 

6 

20 

30 

1 .5446 16388 

4,363323 1 2E  - 02 

6356.813825 

432.4766612 

4.9987043E  - 07 

2.8860033E  - 07 

1,3 

7 

10 

25 

1.527163095 

8.7266463  E - 03 

6356.839812 

648.550316 

3.9224 149E  - 07 

1 .82905 2 1 E - 07 

2.3 

8 

20 

30 

1.527163095 

2.6 1 799388E  02 

6356.839812 

648  550316 

7.496153  E-  07 

4.327906  E - 07 

3,3 

9 

20 

30 

1.527163095 

4.363323 1 3E  - 02 

6356.839812 

648.550316 

7 496153  E - 07 

4.327906  E - 07 

28 


Values  for  the  f"  and  h" 

n m n m 


V 

'o,o 
km  1 

h" 

n0.0 

km  1 

p 

1,0 

km  * 

l»?.o 

km  1 

1 

1.5731 184E-  04 

0 

1 . 5783403 E - 04 

0 

2 

1.5731  1 84E  - 04 

0 

1 .5783403E  - 04 

0 

3 

1.5731  184E  -04 

0 

1 .5783403E  - 04 

0 

4 

1.5731  15  IE  - 04 

0 

1. 5778531 E- 04 

0 

5 

1 .573 1 1 5 IE  - 04 

0 

1 .577853 IE  - 04 

0 

6 

1.5731  15IE-04 

0 

1 .577853 1 E - 04 

0 

7 

1 .573 1087E  - 04 

0 

1.5768788E-04 

0 

8 

1.5731087E  - 04 

0 

1 .5768788E  - 04 

0 

9 

1 .573 1087E  - 04 

0 

1.5768788E-04 

0 

V 

p 

i.i 

km'  1 

hu 
km- 1 

P 

2,0 

km  1 

h2,o 

km  1 

1 

1.3773442E 

-06 

1.2019901  E 

-08 

I.5835193E 

-04 

0 

2 

1.3769246E 

-06 

3.605604  IE 

-08 

1.5835  I93F. 

-04 

0 

3 

1.3760857E 

-06 

6.0081 198E 

-08 

1.5835193E 

-04 

0 

4 

4.1315963E 

-06 

3.6055895E 

-08 

1.5820627E 

-04 

0 

5 

4. 1 3033781' 

-06 

1.08 1 5670E 

-07 

I.5820627E 

-04 

0 

6 

4.127821  IE 

-06 

1.802  24  56E 

-07 

1.58206271 

-04 

0 

7 

6.8845393E 

-06 

6.0080464E 

-08 

1.57915I3E 

-04 

0 

8 

6.8824422E 

-06 

1.8022309E 

-07 

1.5791513E 

-04 

0 

9 

6. 87 8 240 8E 

-06 

3.003 1082E 

-07 

1.5791 5 13E 

-04 

0 

V 

p 

2,1 

km'1 

h2,l 

km  1 

p 

2,2 

km  1 

h2.2 
km  1 

1 

4. 1457488E 

-06 

3.61 79402E 

-08 

3.6I75267E 

-08 

6.3I44163E 

- 10 

2 

4. 14448 5 9E 

-06 

1 .08527 1 8E 

-07 

3.6131  193E 

-08 

1.8935  556E 

-09 

3 

4. 14 1960 7E 

-06 

1.8084 191 E 

-07 

3.6043099E 

-08 

3.I533625E 

-09 

4 

1.24321  20E 

-05 

1.0849347E 

-07 

3.2550933E 

-07 

5.68 1 7864E 

-09 

5 

I.2428333E 

-05 

3.2544735E 

-07 

3.251  1274E 

-07 

1.7038437E 

- 08 

6 

1.2420760E 

-05 

5.42302 10E 

-07 

3.2432006E 

-07 

2.8374329E 

-08 

7 

2.07031 16E 

-05 

1.806733  5 E 

-07 

9.038 143  IE 

-07 

1.5776137E 

-08 

8 

2.0696809E 

-05 

5.4I96503E 

-07 

9.027131 5E 

-07 

4.7309 19  IE 

-08 

9 

2.0684 198E 

-05 

9.030916  IE 

-07 

9.005 1 2 1 7E 

-07 

7.8784606E 

— 08 

Values  for  the  f 


and  h"  (Continued) 


V 

fV 

*3.0 

km'1 

h3,0 

km-1 

"1  1 

E 

h3,l 
km-  1 

1 

1.5886548E-04 

0 

8.3 188628E  - 06 

7.259761 5E -08 

I.5886548E-04 

0 

8.3 163287E  - 06 

2. 1777073E  - 07 

3 

1.5886548E-04 

0 

G.31 12614E-06 

3.628775  IE -07 

4 

1 .5857391 E - 04 

0 

2.4934851  E-  05 

2. 17603 1 5E  - 07 

5 

1 .585739 1 E - 04 

0 

2.4927255E  - 05 

6.5274316E  - 07 

6 

1 .5857391 E - 04 

0 

2.491 2067E  - 05 

1 .0876843E  - 06 

7 

1 .5799 1 54E  - 04 

0 

4. 1485684E  - 05 

3.6204008E  - 07 

8 

1 .5799 1 54E  - 04 

0 

4. 1473047E  - 05 

1.0860100E  -06 

9 

1 .5799 1 54E  - 04 

0 

4. 1 447777E  - 05 

1.8096490E  -06 

V 

f* 

3,2 

knr  1 

h3,2 

knr1 

rv 

3,3 

km  1 

h3.3 
km-  1 

1 

1.8147675E  - 07 

3. 1676885E  - 09 

1 .5834220E  - 09 

4.1463365E-  1 1 

2 

I.8I25565E  -07 

9.4992060E  - 09 

1.579082  E - 09 

1.2427645E-  10 

3 

1.8081372E  - 07 

1 .58191 51 E - 08 

1 .5704 1 38E  - 09 

2.0674889E  - 10 

4 

1 .6324485E  - 06 

2.8494495E  -08 

4.2739029E  - 08 

1.1 191609E  - 09 

5 

1 .6304 596E  - 06 

8.5448767E  - 08 

4.2621885E  -08 

3.3544151E  - 09 

6 

I.6264843E-06 

1.4229893E  - 07 

4.23879 16E- 08 

5.58047  E - 09 

7 

4.529901 8E  - 06 

7.9069730E  -08 

1 97742 13E- 07 

5.1780599E-09 

8 

4. 5 243 82 8E  - 06 

2.371  I286E- 07 

1 .97200I3E  - 07 

1.55  19987E  - 08 

9 

4.5 1 335 1 6E  - 06 

3.94867 10E- 07 

1.961 1762E-  07 

2.5819375E  - 08 

V 

'4,3 

km'  1 

h4 ,3 
km'1 

1 

1.1 120747E-08 

2.91 20701 E - 10 

2 

I.I090266E-08 

8.7282284E  - 10 

3 

1 . 1029387E  - 08 

1.4520463E  -09 

4 

3.0007426E  - 07 

7.857721  IE -09 

5 

2.9925 178E- 07 

2.355 1626E-  08 

6 

2.9760907E  - 07 

3.9 1 80977E  - 08 

7 

1.3875  1 2 2E  -06 

3.6333286E  -08 

8 

1 .3837091 E- 06 

1 .0890027E  - 07 

9 

1.3761  I34E  -06 

1.81 16877E  -07 
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ci  ci  n 


Time-Independent  Coefficients  of  the  Tide  Potential 


n 


0 

I 

I 


3 

3 

3 

3 

4 


m 

“F 

n m 

"F 

n m 

“H 

n m 

"H 

n m 

0 

8.401 8456F.  - 

12 

4.6140106E  — 

12 

0 

0 

0 

8.3681641  E- 

12 

4.5954868E  - 

12 

0 

0 

1 

2.9297877F.  - 

13 

I.6202500E- 

13 

4.3123 569E  - 

15 

2.4544968E  - 

15 

0 

8.3290386E  - 

12 

4. 5 739464 E — 

12 

0 

0 

1 

2.9177584E  - 

13 

1.6I35948E  — 

13 

4.2946458E  - 

15 

2.4444 1 1 8E  — 

15 

2 

2.78405  1 5E- 

15 

I.5423222E- 

15 

8.2 1 32887E  - 

17 

4.683368  E- 

17 

0 

8.2845357E  - 

12 

4.5494  263E- 

12 

0 

0 

1 

2.9046632E  - 

13 

1 .6063487E  - 

13 

4.2753632E  - 

15 

2.4334303E  - 

15 

2 

2.7724443E  - 

15 

1.53589 1 3E  - 

15 

8.1790436E  — 

17 

4.6638388E  — 

17 

3 

1.581308  IE- 

17 

1.0259185E  - 

17 

8. 2068 53  5 E - 

19 

4.6816230E- 

19 

3 

1.8435  105E- 

17 

1.0215973E  - 

17 

8.1722863E- 

19 

4.6619036E  - 

19 

4 -SECOND  COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION 
CAUSED  BY  THE  M2  COMPONENT  OF  THE  OCEAN  TIDE  (BASED  ON  A 
FUNCTIONALIZATION  OF  THE  SEA  SURFACE) 


INPUT  DATA 


C»i 

S„ 


NOAA/Goad  expansion  coefficients  of  the  tidal  surface,  in  meters; 
these  are  specified  by  a computer  card  deck 


R = radius  of  the  earth  (semimajor  axis  of  a suitable  reference  ellipsoid) 
in  kilometers 

suggested  value:  R = 6378.145  km 


G = Newton’s  constant,  in  kilometers,  kilograms,  and  seconds 
suggested  value:  G = 6.6732E-  20  km3  kg  1 sec  2 


/a,  = gravitational  constant  of  the  earth,  in  kilometers  and  seconds 
suggested  value:  pt  = 308601  km3  sec  2 

p = density  of  water,  in  kilograms  and  kilometers 
suggested  value:  p - I.F+  12  kg  km  3 
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o = rate  of  mean  mean  longitude  of  moon 

180 

suggested  value:  o = — 1.405101-  04  deg  sec  1 

7 r 

NMAX  = limit  on  the  number  of  terms  in  the  expansion  of  the  tide  potential; 

NMAX  is  equal  to  the  largest  value  of  the  index  i for  the  expansion 
coefficients  of  the  tidal  surface 

The  parameters  J,  n,  and  t*  specify  the  time  instant  for  which  the  perturbing 
acceleration  due  to  the  tide  is  to  be  found.  Normally,  this  is  the  time  instant  associated 
with  the  current  lime  line  of  the  orbit  integration. 

J = number  of  the  calendar  year 
n = number  of  the  day  within  the  year 
t*  = mean  solar  time  at  Greenwich  (GMT,  UTC),  in  seconds 

X||  inertial,  Cartesian  components  of  satellite  position  vector,  in 

x2  , = kilometers,  associated  with  the  just  specified  time  instant  (for  a 
x3  precise  definition  of  the  reference  frame,  see  Appendices  B and  C of 
Reference  2) 


SUBROUTINE  FOR  THE  TIDE  POTENTIAL  COEFFICIENTS 

This  is  the  first  part  of  the  algorithm  for  the  second  ocean  tide  computer  routine. 
The  time-independent  constituents  of  the  expansion  coefficients  of  the  ocean  tide  potential 
(F”,  F",  H||t  II”)  are  calculated  from  the  expansion  coefficients  for  the  NOAA/Goad  tidal 
ocean  surface.  Once  established,  the  latter  remain  unchanged  as  long  as  the  current  version 
ol  the  NSWC/Schwiderski  Ocean  Tide  Model  (on  which  the  NOAA/Goad  tidal  surface  is 
based)  remains  valid.  Thus,  the  subroutine  under  discussion  will  be  exercised  quite 
infrequently,  namely,  just  once  each  time  a new  edition  of  the  NSWC/Schwiderski  Ocean 
Tide  Model  becomes  available.  That  is  expected  to  occur  at  intervals  of  several  years, 
during  which  the  present  subroutine  will  remain  dormant. 


(4-1) 

(4-2) 

(4-3) 
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(4-4) 


h:; 


2ttR2G p 2 , , 

10  3S  : 


2i+  1 


SUBROUTINE  FOR  THE  PERTURBING  ACCELERATION 


Proceed  exactly  as  in  Equations  3-51  through  3-90  from  the  first  ocean  tide  computer 
routine,  but  delete  Equations  3-53  and  3-54  and  replace  them  with 


Fnm  = Fnm  COS  K + X)  + F"m  Sin  (at*  + X)  (4-53) 

Hnm  = H'nm  COS  K + X)  + H'n'm  sin  (ot*  + X)  (4-54) 

COMPUTER  PROGRAM  OUTPUT 

The  inertial  components  (Tx  ( , TX2,  and  Tx  ^ ) of  the  perturbing  acceleration  due  to 

the  presence  of  the  ocean  tide  are  output  in  terms  of  km  sec-2.  If  required,  convert  to 

m sec-2. 

TRIAL  DATA  FOR  PROGRAM  CHECKOUT 

Non-zero  values  were  assigned  to  selected  NOAA/Goad  expansion  coefficients,  and  the 
corresponding  expansion  coefficients  of  the  tide  potential  were  computed.  The  algorithm  for 
the  perturbing  acceleration  was  then  exercised  with  the  following  results. 


= 1977 

= 202 

= 50000 

sec 

= 6378.145 

km 

= 6.6732E  - 20 

km3 

kg-1 

= 398601 

km3 

sec-  2 
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P = 101  2 kg  km'3 

180 

a - 1 .405 1 9E  - 04  deg  sec" 1 

7 r 

NMAX  = 4 

C20  = +0.2906060089E- 01  m 
C40  = - 0. 107 1 2 1 752E  + 00  m 
C43  = + 0.4357612 19E  - 04  m 

C2'0  = - 0.4424413 130E-  01  m 
C4'0  = + 0.873468034E  - 01  m 
C43  = - 0.160563906E- 02  m 

S20  ' o 

S40  - 0 

s43  = - 0.363303008E  - 02  m 

S20  = 0 
S40  = 0 

S4'3  = - 0.264356490E  - 02  m 

F20  = + 4.9742658E  - 10 
F40  = - 1 .01 86607E  - 09 
F^3  = +4.1438160E-  13 
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I 


F20  = " 7.573211  HE  - 10 
F"0  = + 8.306 13340E-  10 
F"3  = - 1 .526862 IE  - 11 

H'0  = 0 

h;0  = 0 

= - 3. 4547839E-  11 

H"0  = 0 

H40  = 0 

H"3  = - 2.5 1 38645E  - 1 1 

N£  = 933 
AT  = 0.000561 
dQ  = 28325.50056 
T0  = 0.7755099401 
X = 373498.4609 
a = 0.0080511456 
ot*  = 402.557282 
cos(ot*  + X)  = -0.754501 
sin  (at*  + x)  = - 0.656298 

F20  = + 1.2171968E-  10 

h20  = 0 

F40  = + 2.2345060E  - 10 


deg/sec 

deg 
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From  external 


( ABCD)  = 


0 

+ 9.7081 2 16E-  12 
+ 4.2564846F  - 11 

km 
km 
km 

source 


x,  = 3151.52923 
x2  = 5458.60875 
x3  = 3639.07250 


H 


40 


43 


H 


43 


-0.8405285753 
-0.5417605355 
0.26899 13850E-  02 


0.5417623775 

- 0.8405316908 

- 0.5 1 90827376E  - 04 


0.2289080162E  - 02 
0.1413662999E-  02 
0.9999963803 


y,  = + 3 16.6486 1 km 

y2  = - 6290.36338  km 

y3  = + 3647.25332  km 


3^20 

3y, 

3V2q 

3y, 

3U20 

3y2 

**20 

3y2 


= - 0.0000964047 


= - 4.7E  - 13 


= +0.0019151219 


= - 2.4E  - 14 


km/sec2 


kni/sec2 


km/sec2 


km/sec2 
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Values  for  the  Eigenfunctions  of  the  Tide 
Potential,  in  km2  sec"2 


n 


m 


U 


n m 


o -1 

0 0 


0 

54.76683963 


1 

1 

1 


-1 

0 

1 


-1.044042677 
24.051 191 13 
2.088085354 


-1  -0.4584976993 

0 -5.186455141 

1 2.750986196 

2 -94.01442066 


-1 

0 

1 

2 

3 


-0.0512402699 

-16.10992273 

0.6148832385 

-206.435027 

-53.85812191 


-1 

0 

1 

2 

3 

4 


0.109342534 

-9.393545781 

-2.18685068 

-136.7982658 

-165.56486 

1863.67849 


-1 

0 

1 

2 

3 

4 


0.0917032956 

2.472201758 

-2.751098868 

136.846714 

-182.4236505 

7366.011828 


-1 

0 

1 

2 

3 

4 


0.0153004451 
8.002095405 
-0.6426186943 
349.1 179421 
45.32033435 
11350.89177 


V 

n m 


o 

o 

-20.74036524 

0 

-41.48073047 

-9.108257692 

0 

-54.64954615 

-9.489169689 

-1.017910413 

0 

-12.21492495 

-20.83613333 

354.226451 

2.172137348 

0 

43.44274695 

-13.80747707 

1088.924921 

380.085929 

1.821726148 

0 

54.65178443 

13.8123671 
1 199.805712 
1502.253452 

0.303950046 

0 

12.76590193 

35.23756643 

-298.0731712 

2314.94556 
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™2  0 

9y3 

av20 

9y3 

au4o 

9y, 


9V4q 

9Vi 

9U40 

9y2 

9V4  0 

ay2 


9U4q 

9y3 

9V4q 

9y3 

9U43 

ay. 


a v 


4 3 


9y, 


9U43 

9y, 


9V43 

9y2 


= +0.0075774019 

= 0 

= +0.0004313321 

= +7.8E-  13 

= - 0.0085686018 

= 0 

= - 0.0019380257 

= 0 

= - 0.5130748473 

= - 0.1 1 12689700 

= - 0.1242624348 

= + 0.64 1808  246  i 


km /sec2 


km /sec 2 


km /sec 2 


km /sec 2 


km /sec 2 


km /sec2 


km/sec2 


km /sec 2 


km /sec2 
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3U4  3 

= +0.0572027292  km/sec2 

3y3 

3 V4  3 

=-  0.3762240313  km/sec2 

0y3 

Ty  = - 9.632495 E-  12  km/sec2 

Ty2  = + 2.443056E  - 11  km/s^c2 

Tyj  = - 1.496932  IE  - 11  km/sec2 

Tx  = - 5.179392E  - 12  km/sec2 

TXj  = - 2.5752406E-  11  km/sec2 

TXj  = -1.495676E-  11  km/sec2 


5 - COMPUTER  ROUTINE  FOR  THE  PERTURBING  ACCELERATION 
CAUSED  BY  THE  TIDES  OF  THE  SOLID  EARTH 


INPUT  DATA 

The  parameters  J,  n,  and  t*  specify  the  time  instant  t for  which  the  perturbing 
acceleration  due  to  the  tide  is  to  be  found.  Normally,  this  is  the  time  instant  associated 
with  the  current  time  line  of  the  orbit  integration. 

J = number  of  the  calendar  year 
n = number  of  the  day  within  the  year 
t*  = mean  solar  time  at  Greenwich  (GMT,  UTC),  in  seconds 

At  = tidal  lag  parameter,  in  seconds 

w = earth’s  sidereal  rate  of  rotation 

suggested  value:  to  = 4.178074622E  - 03  deg  sec- 1 


39 


K 


l 

i 


R 


e 
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the  inertial,  Cartesian  coordinates  of  the  moon,  in  kilometers,  for  the 
= time  instant  (t-  At),  from  the  TERRA/CELEST  moon  ephemeris;  the 
reference  trame  is  relerred  to  mean  equator  and  mean  equinox  at  time 
instant  t, 

the  inertial,  Cartesian  coordinates  of  the  sun,  in  kilometers,  for  the  time 
= instant  (t  - At),  trom  the  IERRA/CELEST  sun  ephemeris;  the  reference 
frame  is  relerred  to  mean  equator  and  mean  equinox  at  time  instant  tj 

= epoch  time  for  the  moon  ephemeris  and  sun  ephemeris,  in  Julian  Date 
t,  = JD2433282.423 

inertial.  Cartesian  components  of  satellite  position  vector,  in  kilometers, 
= associated  with  the  time  instant  t(J,  n,  t*)  (for  a precise  definition  of 
the  reference  frame,  see  Appendices  B and  C of  Reference  2) 


_ Cartesian  satellite  velocity  components  associated  with  time  instant  t,  in 
km  sec"1 

= gravitational  constant  ot  the  earth,  in  kilometers  and  seconds 
suggested  value:  pE  = 398601  km3  sec"2 

= radius  of  the  earth  (semimajor  axis  of  a suitable  reference  ellipsoid)  in 
kilometers 

suggested  value:  R = 6378.145  km 

= square  of  eccentricity  of  reference  ellipsoid 
suggested  value:  e2  = .006693421623 

in  case  it  is  desired  to  start  from  the  flattening,  f,  find  e2  from 
e2  = (2  - Of 

= expansion  coefficients  for  the  Love  numbers 
suggested  values:  k20  = 0.3 
k2,  = 0.01 
k3Q  - 0.1 
k3|  = 0.0 1 
k22  = 0.1 
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m'  = mass  of  moon  or  sun,  in  kilograms 

suggested  values:  m'L  = mass  of  moon  = 7.369328 IE  + 22  kg 
m^  = mass  of  sun  = 1.99E  + 30  kg 

M = mass  of  earth,  in  kilograms 

suggested  value:  M = 5.973 16 13E  + 24  kg 

a'  = semimajor  axis  of  moon  or  sun  orbit  relative  to  earth,  in  kilometers; 
select  values  that  are  compatible  with  the  moon  ephemeris  or  sun 
ephemeris  (optionally,  values  may  be  obtained  from  a current  edition  of 
the  AMERICAN  EPHEMERIS) 
suggested  values:  a'L  = 384400  km 

a's  = 149500000  km 


SUBROUTINE  FOR  THE  POSITION  OF  THE  FICTITIOUS  MOON 


(5-1) 


Xx  = cos  £0  cos  0 cos  z - sin  £0  sin  z (5-2) 

Yx  = - sin  £0  cos  0 cos  z - cos  £0  sin  z (5-3) 

Z - - sin  0 cos  z (5-4) 

Xy  = cos£0  cos  0 sin  z + sin  cos  z (5-5) 

Yy  = - sin  £0  cos  6 sin  z + cos  £0  cos  z (5-6) 

Z - - sin  0 sin  z (5-7) 

Xz  = cos  £0  sin  6 (5-8) 

Yz  = - sin  £0  sin  0 (5-9) 

Zz  = cos  0 (5-10) 
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(n  = 1(2304.250  + 1.396Tn)T  + 0.302T2  +0.018T3(  — — 
0 ' 0 ’ 3600 

z = + j 0.79  IT2  | 

0 ' ’ 3600 

6 = ] (2004.682  - 0.853Tn  )T  - 0.426T2  - 0.042T3  ( — — 

> u > iccnn 


(5-11) 


(5-12) 


(5-13) 


£0,  z,  and  0 result  in  degrees. 


T0  = (t,  - 2415020. 313)(0.273790926497E  - 04) 


(5-14) 


T = (0.273790926497E  - 04)(t2  - t,  ) 
(t2  - t,  ) = ) (J  - 50)365  + A0  + (n  - 1)  + .077  \ 


(5-15) 


(5-16) 


/J  - 53 

A,.  = INI  + 1 


(5-17) 


The  Cartesian  coordinates  of  the  fictitious  moon  (x*,  y*,  and  z*)  referred  to  the 
coordinate  system  in  which  the  satellite  orbit  integration  is  performed  may  rtow  be 
computed: 


cos  w At 
sin  co At 
0 


-sin  <o At  0 


cos  coAt  0 II  yj^ 
0 \J\z[ 


(5-18) 


Also,  compute 


= + v/(x*)2  + (y*  )2  + (z*  )2 


(5-19) 


(5-20) 


(5-21) 


z* 

v.._c 


(5-22) 
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SUBROUTINE  FOR  THE  POSITION  OF  THE  FICTITIOUS  SUN 


(5-23) 


Evaluate  the  matrix  elements  Xx,  Yx>  etc.,  according  to  Equations  5-2  through  5-17. 
Then,  calculate  the  Cartesian  coordinates  of  the  fictitious  sun  (x|,  y*,  z*)  referred  to  the 
coordinate  system  in  which  the  satellite  orbit  integration  is  performed: 


cos  coAt 
sin  cuAt 
0 


-sin  cuAt 
coscuAt 
0 


(5-24) 


Further,  compute 


+ V(X*)2  +(y*)2  + (z*)2 
X? 


SUBROUTINE  FOR  SATELLITE  POSITION 


(5-25) 

(5-26) 


(5-27) 


(5-28) 


r = +Vx?  +x$  + x} 


(5-29) 


X = 


*l 

r 


(5-30) 
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M = 

^2. 

(5-31) 

v = 

r 

*3 

(5-32) 

ii 

tN 

> 

r 

(x,  )2  + (x,)2  + (x3)2 

(5-33) 

a = 

(--  -)" 

(5-34) 

\ r Pf  / 

SUBROUTINE  FOR  THE  PERTURBING  ACCELERATION 


Tx  t , TX7,  and  TX;j  are  the  Cartesian  components  of  the  perturbing  acceleration,  in 
km  sec- 2.  These  are  referred  to  the  (inertial)  reference  frame  in  which  the  satellite  orbit 
integration  is  performed,  and  are  ,hus  readily  suited  for  insertion  into  the  TERRA/CELEST 
equations  of  motion. 


Tx,  = 7]  _Co  + + 7T(5V2-2p;)+  7(7V3-F)+  ^(9V4-H)]> 


c,  , c2  c3  c4 

~a3  + 7^11  + 7^12  + 7 p 


(5-35) 


x 2 


7ilco  + 


3C,V,  C 


j-1  + rT(5V2-2p;i)+  7^3  -F)  + 7(9V4-H) 


Cl  , C2 

t-a;  + — p 


2 I 


C3  c4  ] ( 

r2  ^2  2 + r3  p23  J ( 


■ ?![c»+ 

f-  4C,  , C2  C3  C4  1 I 

' 5 A0  + r ^3  1 + r2  ^3  2 + r3  ^3  3 | 


(5-36) 


3-£71  + 7"  (5V2  - 2Pq)  + ^f(7V3-  F)+  7(9V4-  H) 


■J 


(5-37) 


(5-38) 


pM  = 2\p;  +mp;  + < 

p,  2 = ( i - 5«/J  >s;  + 3(\2  - m2  >s;  + 2X(3ais;  + «^s;  > + ^s;  <s-39) 

p, j = v (\  - ^2) r;  + 3wx2  - mj >t;  + 6xm«/t;  + ( i - iv1  x2xt;  + mt; » <5-40  > 

p2l  = -2pp;  +xp;  +<  (5-4D 

p22  = (i  - 5»»j)s;  - ftx^s;  + 3<x2  - m2)s;  +^xs;  - 2^)  (5-42) 

p2  3 = »-/»-  j r;  - 6x^t;  + 3^x2  - m2  it;  + < 1 - iv2  ><xt;  - :^t'6  > (5-43 ) 

Pj,  = - (h;p;  +xp;  +#ip;  (5-44) 

(>J2  = - iowxs;  +jus2)  + (3- is^is;  +<x2  - fi2)s'b  +xps2  (5-45) 

pj  J = ( 1 - Iv1  ((XT',  + uT2  ) - 20i43  - lv2  )T; 

+ X(X2  - 3/u2  IT;  + ^(3X2  - H*  )T'S  (5-46) 

i4H<x2  - m2  )t:  +x#jt;i 

1 = 2(xs;  +ps;  + 3^s; » (5-47) 

h = 2 [hxt;  + n \'2 ) + 6<  1 - 5i>2  »r,  + (x2  - n2  >t;  + xpt;]  (5-48) 


where 


Co 

^3  E ^2  0 <■  ^2  2^  UKA0 

(5-49) 

= Kak2  ( J2 

(5-50) 

r2 

= Ka2aJ 

(5-51) 

= KaJa4 

(5-52) 

(4 

= Ka4as 

(5-53) 
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and 


and 


m 


k = — n2  it2a'30*3 
M 


and 


a3  n2 


fj  = M(I 
R 


a = — 


« = — 


a*  = — 


Vo  s 


v,  = ;(A'»X  + A>  4A<>» 


V,  = 


4 

Z !>'  I> 

..  II  II 


V,  = 


Z S'S 
1 " " 


V,  = 


Z T'  T 

n I " n 


A'  = 


- ( I - h**1 ) 
4 


A'  = 


- (X*2  - #j*2  ) 
4 


A'  = 


A',  = 


A'  = 


3X*/i* 

3**v* 

3fj*v* 


(5-54) 

(5-55) 

(5-56) 

(5-57) 

(5-58) 

(5-59) 

(5-60) 

(5-61 ) 

(5-62) 

(5-63) 

(5-64) 

(5-65) 

(5-66) 

(5-67) 
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▲ 


(5-68) 


and 


- A*(  I - 5u*2 ) 
8 


(5-69) 


~/u*(l  - 5u*2 ) 

O 


(5-70) 


- u*(3  - 5u*2 ) 
4 


(5-71) 


- X*(X*2  - 3m*2  ) 
8 


(5-72) 


- m *( 3X *2  - n*2 ) 
8 


(5-73) 


u*( A*2  - /i*2 ) 


(5-74) 


= 15A*/a*n* 


(5-75) 


- (k,o + |k=»)(i  - if!)A;  + 


f.  - I - 3-' 


K • (k;o  ' 7ki!)(|-  ^’)A',  * 7*3,“'9*BL 


P,  = A2  - M2 


n - (k.°-  fk»)(>-  i7,!)A;  + vk..“,‘j,B; 


(5-76) 


(5-77) 


(5-78) 


(5-79) 


(5-80) 


P2  = Am 


^ * (k’°t7Lk’>)(|-  ii*0A»-  fk>.",(J*B; 


(5-81) 


(5-82) 


P3  = Av 


■ (k,.  + f‘,1)('-  He!)A:- 


(5-83) 


(5-84) 


P4  = nu 


(5-85) 
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a nil 


S',  ="  jl41A;+kJ#a'P*B'|  (5-86) 

s,  = Ml  - 5^2)  (5-87) 

s;  =- jk21A;+k30a'^B;  (5-88) 

s2  = n(\  - 5v2 ) (5-89) 

S3  = | k2 . A0  +k30CV'^B3  (5-90) 

53  = v(3  - 5v2 ) (5-91) 

54  = k3oa'^*B4  (5-92) 

S4  = X(X2  - 3ju2 ) (5-93) 

s;  = k30^*B_;  (5.94) 

Ss  = /a(3X2  - /a2 ) (5-95) 

= k2i  A'.  +k3oa'(3*Bl  (5-96) 

56  = n(X2  - n2 ) (5-97) 

57  = k 2 1 A 2 + k3oa'^B7  (5-98) 

S7  = Xw  (5-99) 
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and 


n ={  7ik”)'  r4k”}A«  + 

t,  -w(i- 1*1) 


1 3 

/ 2 \ 

9 \ 

k + — k , ) e2  - 

— k 

1 !4 

\ 2 O -j  22  1 

70  22 

T3  = 3 - 30u2  + 35y4 

T;  = k31a'/3*B; 

T4  = MX2  - 3/a2) 

t;  = k3I«^*B; 

T$  = nv(3\2  - M2 ) 

T‘  " (fi(k»o- 

r6  = (X2  - m2ki  - 7t-2) 

n * {iHk=°-  7k”) 

T,  =MI-  7t>2) 


(5-100) 

(5-101) 

(50102) 

(5-103) 

(5-104) 

(5-105) 
(5-106) 
(5-107) 
(5-108) 
(5-109) 

V,  - ^k31a'/?*B;  (5-110) 

(5-111) 

A2  ~ 7k3.tt'0*B7  (5-112) 

(5-113) 


Note  that  the  components  Tx  of  the  perturbing  acceleration  depend,  essentially,  on 
the  parameters 


X*,  /a*.  v *,  r*.  a',  m',  x,,  x2,  x3  , 
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Tv  = Tx  (X*,  p*,  i>*,  r*,  a',  m',  x,,  x,,  x, 

* i * i 12  3 


To  obtain  the  perturbing  acceleration  due  to  the  moon  tide,  substitute  for  these  parameters 
those  specifying  lunar  position,  orbital  semimajor  axis,  and  mass. 


(TXj)L  = Tx  .(X* , m*, 


rL  ’ aL  ’ 


m. 


X3) 


Otherwise,  substitute  the  solar  parameters  to  produce  the  perturbing  acceleration  caused  by 
the  tidal  action  of  the  sun. 

(Tx.)s  = Tx.(X*,  f's’  rs ’ as 1 ms ’ Xl’  x2 ’ X3 ^ 

COMPUTER  PROGRAM  OUTPUT 

The  inertial  components  (TX()L,  (TX,)L,  (TX3)L,  (TX])S,  (TXi)$,  and  (TX3>S  of  the 
perturbing  acceleration  due  to  lunar  and  solar  tides  of  the  solid  earth  are  input  in  terms  of 
km  sec"2. 


TRIAL  DATA  FOR  PROGRAM  CHECKOUT 

The  following  trial  case  is  restricted  to  the  lunar  tide.  While  executing  this  case,  the 
entire  algorithm  was  exercised  and  it  was,  thus,  unnecessary  to  also  consider  the  solar  tide. 


Input  Data 


t •••  J = 1977 

n = 88 

t*  = 57000 

sec- 

At  = 100 

sec 

w = 4.178074622E-  03 

deg  sec 
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Epoch  for  the  moon  ephemeris  (reference  frame  2j) 
t,  = 1950.0  = JD2433282.423 
(1950.0  means  beginning  of  Besselian  year  1950) 

Epoch  for  the  satellite  orbit  (reference  frame  2,) 


t2-”  J = 1977 
n = 88 

t*  = 0 sec 

t2  = JD244323 1.5 

Moon  coordinates,  in  frame  at  time  (t  - At) 
x"  = - 184258.5823  km 

Yl  = 329791.1351  km 

z"  = 103840.2349  km 

Satellite  position  and  velocity,  in  frame  Z2  , at  time  t 
x,  = - 4009.582237  km 

x2  = + 103.9008135  km 

Xj  = - 5269.570696  km 

Xj  = - 5.978210289  km  sec'1 

x2  = + 1.710442216  km  sec-1 

x3  = + 4.584971066  km  sec'1 

= 398601  km3  sec'2 

R = 6378.145  km 


e2  = 6.69342 1623E  - 03 

k20  - 0.3 
k2,  = 0.01 
k30  = 0.1 
k3]  = 0.01 
k22  = 0.1 

m'L  = 7.3693281  E + 22  kg 

M = 5.973 1 6 1 3 E + 24  kg 

a'L  = 384400  km 
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Results 


0 = 0.1516444800 

deg 

£ = 0.1744119454 

deg 

z = 0.1744282487 

deg 

'0  = 0.5000000017 

centuries 

T = 0.2723967010 

centuries 

Xx  = + 0.9999779632 
Yx  = - 0.0060883617 
Zx  = - 0.0026466801 
Xy  = +0.0060883617 
Y = +0.9999814657 
Zy  = - 0.0000080574 
Xy  = +0.0026466801 
Yz  = - 0.0000080567 
Zz  = +0.9999964975 


\'L  = - 186537.2414 

km 

y'L  = +328662.3531 

km 

z^  = + 103349.5407 

km 

r*  = 391785.9267 

km 
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x*  = - 188928.9046 


km 


y*  = + 327293.3757 
z*  = + 103349.5407 

X*  = - 0.4822248369 
/i*  = + 0.8353882910 
p*  = +0.2637908451 
X = - 0.6054593779 
M = +0.0156893457 
v = - 0.7957215507 


r = 6622.380268 

km 

a = 6567.451 187 

km 

K = 0.323073882E  - 05 

km2 

n = 0.1 18624368E  - 02 

sec' 

a'  = 0. 1659246878E  - 01 
a = 0.971 1750904 
li*  = 0.981 1480551 

A„  = +0.1978107925 
A'  = - 0.3489996026 
\'2  = - 0. 1 208534947E  + 0 1 
A 3 = - 0.3816194918 
a;  = +0.661  1033498 
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km 

km 


B',  = - 0.1  179169837 
B'2  = +0.2042749770 
B'3  = +0.1748980753 
B;  = + 0.56091  18737 
Bj  = - 0.000131 1647178 
Bg  = - 0.4603145005 
B2  = - 0. 1 594002275E  + 01 

Po  = +0.06443748370 
P0  = - 0.5995183587 
P',  = - 0.09451271981 
P,  = +0.3663349027 
P'2  = - 0.3272838250 
P2  = - 0.009499261487 
P'3  = - 0.1 190554808 
P3  = + 0.4817770751 
P;  = +0.2062472668 
P4  = - 0.01248435049 

S'  = +0.0005712740433 
S,  = +0.131  1 342628H  + 01 
S',  = - 0.0009896538092 
S2  = - 0.03398098796 
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s'3  = +0.001471593023 

53  = +0.1319815046 

S;  = +0.0009131459348 

54  = - 0.2215028279 

Sj  = - 0.2 1 353 1 8 1 2 1 E - 06 

55  = + 0.0172503888 

Sj.  = - 0.004239372772 

56  = - 0.2915005769 
S'7  = - 0.01468033233 

57  = + 0.007558767081 

T'  = +0.02363141 132 
T,  = - 0.2300019019 
T'2  = - 0.04093817407 
T,  = +0.005960068473 
T'j  = - 0.002454126571 
T3  = - 0.196341 1384E  + 01 
T;  = + 0.9 13  I459348E  - 04 
T4  = + 0.1762545737 
T'j  = - 0.21 353  18 1 2 1 E - 07 
Ts  = - 0.01372650615 
= +0.02595922645 
T6  = - 0.1 257338 135E  + 01 
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N 


T’7  = +0.08989303178 
T?  = +0.03260345555 

Vn  = + 1 .0000000000 
V,  = +0.1742073243 
V2  = - 0.1494101 174 
v3  = +0.1 8995 34786E-  02 
V4  = - 0.03055341082 


= - 0.0001622651725 

km3 

_ 2 
sec 

= + 0.135329691 5E  + 01 

km4 

sec  * 

= + 0.863  1523953E  + 06 

km5 

- 2 
sec 

= + 0.550531 1134E+  10 

km6 

- 2 

sec 

= +0.351  1367268E+  14 

km7 

2 

sec  * 

F = - 0.007748690189 
H = +0.1048876747 


p,  , = + 0.2040473678 
p, 2 = - 0.004135328808 
p,  3 = + 0.1119464221 
p2|  = +0.03700735158 
p22  = - 0.0049832331 13 
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p23  = + 0.1740434693 
p31  = +0.3829649087 
p32  = - 0.01385120787 
p33  = + 0.2036553308 


- 0.1929747594E  - 

09 

km 

sec  2 

+ 0.9604765343E  - 

10 

km 

sec-2 

- 0. 1 824727285E  - 

09 

km 

sec"2 

T = 0.2824 193799E  - 09  km  sec"2 


DEVELOPMENT  STATUS  OF  THE  COMPUTER  ROUTINES 


By  early  summer  1978,  all  five  tide  routines  had  been  coded  and  made  part,  as 
operational  subroutines,  of  the  TERRA  and  CELEST  computer  program  systems.  Program 
checkout  was  accomplished  (separately  for  each  tide  routine)  by  first  selecting  typical  input 
data  that  would  exercise  the  algorithm  to  a sufficient  degree.  Subsequently,  the  programmer 
would  run  the  algorithm  under  test  on  the  computer,  while  the  author  would  perform 
identical  calculations  on  a programmable  electronic  calculator.  Agreement  to  ten  significant 
digits  of  the  computer  results  with  those  of  the  manual  computation  was  taken  as  evidence 
for  successful  computer  implementation  of  the  algorithm. 

As  indicated  by  the  last  paragraph  of  the  INTRODUCTION,  the  two  ocean  tide 
algorithms  may  be  expected  to  produce  nearly  identical  perturbing  accelerations  because  the 
tidal  ocean  surface  of  the  second  algorithm  is  a functionalization  of  the  discrete 
representation  of  the  ocean  surface  on  which  the  first  algorithm  is  based.  Like  the  routines 
for  the  two  air  tides  and  that  for  the  solid  earth  tide,  each  of  the  two  ocean  tides  had 
been  checked  out  individually.  In  the  process,  to  facilitate  the  manual  calculations,  a simple 
ocean  surface  had  to  be  selected  for  each  of  the  two  algorithms.  Because  of  the  different 
mathematical  structure  of  the  two  algorithms,  simplicity  meant  that  each  algorithm  required 
its  own  individual  surface.  Thus,  it  was  impossible  to  directly  compare  the  resulting 
perturbing  accelerations.  This  in  turn  meant  that,  during  this  phase  of  the  checkout,  the 
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mutual  consistency  of  the  two  ocean  tide  routines  could  not  be  verified.  Because  the  two 
routines  are  so  complex,  this  remained  a source  of  concern  despite  the  caution  practiced 
while  deriving  them.  The  question  was  finally  settled  by  running  each  algorithm,  using  the 
complete  input  data  base.  The  resulting  perturbing  accelerations  agreed  with  each  other  to 
20  percent.  Considering  the  differences  in  the  physical  definitions  of  the  tidal  mass  layers 
(one  being  a discrete  point  set,  and  the  other  being  a continuous  surface),  this  was 
interpreted  as  satisfactory  evidence  for  the  validity  of  the  algorithms. 

For  operational  use.  that  ocean  tide  routine  was  selected  which  employs  the  discrete 
tide  surface  representation.  Realizing  that  the  second  algorithm  contains  a more  accurate 
physical  model  (by  virtue  of  its  continuous  tidal  mass  layer  as  compared  with  the  point 
mass  approximation  inherent  in  the  first  ocean  tide  routine),  it  was  yet  considered  a 
decisive  advantage  of  the  algorithm  chosen  that  it  is  capable  of  accepting  among  its  input, 
and  of  being  capable  to  directly  process,  without  the  need  for  any  intervening 
computational  steps,  the  earlier  mentioned  tide  amplitude  and  phase  tables  that  are  in-house 
products  of  NSWC.  This  aspect  is  particularly  important  because  it  is  expected  that,  during 
the  useful  life  ot  TFRRA/CELEST,  the  present  M,  tide  tables  will  be  augmented  by 
similar  tables  for  various  other  constituents  of  the  ocean  tide  spectrum.  Once  available,  a 
mere  extension  of  the  present  algorithm  will  make  possible  the  immediate  use  of  the 
additional  program  input,  without  first  having  to  obtain  the  additional  functionalizations  of 
the  tidal  ocean  that  would  be  needed  if  the  second  ocean  tide  routine  had  been  adopted 
for  inclusion  in  TERRA/CFLEST. 
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NOTES  CONCERNING  THE  SECOND  OCEAN  TIDE  ROUTINE 


As  mentioned  in  the  INTRODUCTION,  the  second  computer  routine  for  the  ocean 
tide  is  based  on  a functionalization  of  the  tidal  ocean  surface.  This  functionalization 
represents  the  height  of  the  tidal  surface  above  mean  sea  level,  h,  as  an  expansion  in 
surface  harmonics. 


h = f cos  (ot*  + x ~ 5)  = ? cos  6 cos  (a t*  + x)  + f sin  6 sin  (at*  + x)  (A-l) 


1 5 n 

f cos  6 = £ £ (Cnm  cosmX  + Snm  sin  mA)P'nn  (sin  0) 

n = 0 m =0 


1 5 n 

f sin  5 = 52  5"  (C’  cos  mX  + S'  sin  mX)Pm  (sin  0) 

L-*  /L~i  nm  nm  n v ' 

n = 0 ni  = 0 


where 

X = longitude 
0 = latitude 

C , S„m,  C'  S'  = expansion  coefficients 

nm  nm  nm  nm  ^ 

a = rate  of  the  mean  mean  longitude  of  the  moon 
X = mean  mean  longitude  of  the  moon,  at  the  beginning 
of  the  day  UT 

Note  that  h is  essentially  a function  of  longitude,  latitude,  and  time.  This  function  is  valid 
everywhere  on  the  globe  (sea  and  land  - see  the  comments  at  the  end  of  this  appendix). 

The  surface  density  corresponding  to  the  tidal  mass  layer  is 


1 ( | 
s = ph  = — -|0  cos  (ot*  + x)  + (3  sin  (at*  + x)^ 


a = Gp  f cos  6 


(3  = Gp  f sin  8 
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The  gravitational  potential  associated  with  the  tidal  bulge  is  then 


s dS 

<p  = Gff— — — . = 0,  cos  (at*  + x)  + 0,  sin  (at*  + x) 
|r  - r | ' 1 
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Now,  evaluate  0. 


J2  n ..  + n/2  / 

dX'  ( cosO'dd'  ' 

o J-„/2  /IF- r 


(A-12) 


Let 


sin  0 = y ( A- 1 3 ) 

sin  O'  = y'  (A-14) 


ami  consider  the  term  ij  from  the  expansion  for  f cos  6 , namely 


(f  cos  6)(j  = (C4j  cos  jX'  + Stj  sin  jX')pj(y') 


Then, 
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(A- 16) 


Repeated  application  of  the  orthogonality  relationships  for  the  Legendre  functions,  as 
well  as  for  the  sine  and  cosine  functions,  makes  possible  the  following  equations: 


>2"I 


+ 1 


(0i  ,.J=o  " GpR2Cjo27r  j dy'  P,(y')  ~ Pn(y)  p (/) 

n = 0 r 

R‘ 


= 2»r2gp  2-^  -r,  c,„  P,<y) 


(A- 17) 
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sin  jX'  cosm(X-  X')  dX'J 

J2ir 

cos  jX  cos  m(X  ~ X)  dX;  = 7r  5^.  cosjX 


f2- 
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sin  jX'  cosm(X-X')  dX'  = n 8 . sin  iX 

mj  J 


1 n = m 

nm  I ii* 

0 n # m 

2 R* 


(0i  \j*  o 2^pR2  2—  P](sin  0)  * (Cy  cosjX  + Sy  sin  jX) 


and  consequently, 


2 Ri 

</>,  )„  - 2ttR2Gp  2~  — j * (Cy  cosjX  + Sy  sin  jX)  Pj(sin  0) 


Thus,  the  following  results  in  a similar  fashion: 
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2ttR2Gp  2 r , 

IS  cos(ot*  + x)  + S sin  (at*  + x) 

^ 2v  + 1 *-  "M  *'M  -l 

R1’ 

% rTTT  (sin  0)  cos  pX 


(A-27) 

(A-28) 


V.-4  = ?v  (sin0>  sin^  (A-29) 

Up  to  here,  each  of  the  physical  quantities  occurring  in  the  equations  may  still  be 
entered  in  terms  of  arbitrary  units.  However,  for  the  purpose  of  the  second  ocean  tide 
computer  algorithm,  it  is  necessary  to  employ  compatible  units:  kilometers  for  length, 
kilograms  for  mass,  and  mean  solar  seconds  for  time.  As  obtained  from  NOAA,  the 
coefficients  C^,  C'  , and  are  listed  in  meters.  To  accommodate  the  existing 

computer  card  deck  for  these  data,  it  is  necessary  to  multiply  the  right  sides  of  Equations 

A-26  and  A-27  by  10~3  prior  to  their  use  in  the  main  part  of  this  report. 

It  should  finally  be  mentioned  that  the  expansion  for  the  tidal  surface  (Equations  A-l 

through  A-3)  is  intended  to  be  valid  for  both  sea  and  land  in  the  sense  that  it  is  a fit  to 

the  NSWC/Schwiderski  tidal  heights  and  phases  on  the  ocean  and  to  zero  tidal  heights  on 
land.  Because  of  the  relatively  limited  number  of  terms  employed  in  the  expansion,  this  fit 
is  somewhat  imperfect,  especially  near  the  coast  lines.  This  imperfection  is  intentional.  The 
reason  for  it  is  that  the  expansion  was  originally  meant  to  be  the  basis  for  finding  a tide 
potential  for  satellite  work  rather  than  a continuous,  highly  accurate  description  of  the 
tidal  surface.  Satellite  orbits  appear  to  be  much  more  sensitive  to  the  low-order  harmonics 
in  the  tidal  disturbing  potential  than  to  the  higher  ones.  Based  on  information  from  an 
outside  source,  we  in  fact  anticipate  that  certain  of  the  second-  and  fourth-order  harmonics 
will  be  the  main  contributors  to  the  orbital  perturbation.  Thus,  the  emphasis  on  the 
low-order  terms  apparent  in  the  expansion  for  the  tidal  surface  is  rather  an  advantage  than 
a defect  as  far  as  our  own  satellite  work  is  concerned. 
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A COMPARISON  OF  THE  PERTURBING  ACCELERATION  DUE  TO  THE 
SOLID  EARTH  TIDE  WITH  THAT  CAUSED  BY  THE  OCEAN  TIDE 


Preliminary  runs  of  the  first  computer  routine  for  the  ocean  tide  and  of  the  routine 
tor  the  solid  earth  tide,  all  of  which  used  complete  sets  of  input  data  (as  opposed  to  the 
restricted  input  data  of  the  numerical  trial  cases  listed  in  the  main  body  of  this  report), 
have  resulted  in  the  following  perturbing  accelerations. 

The  satellite  coordinates  and  the  perturbing  acceleration  components  listed  are  referred 
to  the  inertial  frame  corresponding  to  the  mean  equinox  and  mean  equator  of  J = 1977, 
n = 202,  and  t*  = 0 sec. 


Time  t 


J = 1977 
n = 202 
t*  = 50000  sec 

Satellite  coordinates  at  t 

x,  = 3151.52923  km 
x2  = 5458.60875 
x3  = 3639.07250 


Perturbing  acceleration  due  to  ocean  tide  at  t 


Tx  - + 0.594629E  - 11  km  sec  2 
Ty  = + 0. 134622E  - 10 
Tz  = - 0.258193E-  11 


Perturbing  acceleration  due  to  solid  earth  tide  at  t 


Tx  = + 0. 149368E  - 09  km  sec~2 
Ty  = + 0. 125475E  - 10 
Tz  = - 0. 184573E  - 10 

Tx  = - 0.493523E-  10 
Ty  = +0.38I372E-  10 
Tz  * + 0.713397E  - 12 


1 V 

> due  to  moon 

' ( Nominal  values  for 

ILove  numbers 

due  to  sun  > 
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Tx  = + 0.1 64304 H - 09 
Ty  = + 0.013802E  - 09 
Tz  = -0.020303E  - 09 

Tv  = -0.054288E  - 09 
Ty  = + 0.041951  E-  09 
Tz  = + 0.000785  E-  09 

Tx  = + 0.1 7924 IE  - 09 
Ty  = + 0. 150570E  - 10 
Tz  = - 0.22 1487E  - 10 

Tx  = -0.592228E-  10 
Ty  = + 0.457647E-  10 
Tz  = + 0.856076E  - 12 


The  magnitudes  of  the  perturbing  accelerations  are  then 

TfOcean  tide)  = ^ocean  = 1-49E-  11  km  sec-2 

TfSolid  earth  tide  due  to  moon,  for  nominal  Love  Numbers) 
= TSOL|D  = 1.51026E-  10  km  sec”2 

T(Solid  earth  tide  due  to  moon,  for  Love  numbers  increased 
= TSOLiL)  + )0  = 1.66048E-  10  km  sec”2 

TfSolid  earth  tide  due  to  moon,  for  Love  numbers  increased 
= Ts o L i d + 2 o = L81231E  10  km  sec’2 

Note  that 


due  to  moon 


due  to  sun  / 


due  to  moon 


due  to  sun 


or 


^SOLID+IO  OLID  31  ^OCEAN 


^SOLlD  + 10  ~ ^SOLID  ^ ' I^OCEAnI 
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Each  Love  number  , 

increased  by  10  percent 

i 


/ 

Each  Love  number 
increased  by  20  percent 


by  10  percent) 

by  20  percent) 


(B-l ) 


(B-2) 


t 


This  indicates  that,  for  the  particular  combination  of  time  and  satellite  position 
chosen,  if  one  increases  the  Love  numbers  by  ten  percent,  the  perturbing  acceleration  due 
to  the  lunar  component  of  the  solid  earth  tide  increases  by  an  amount  roughly  equal  to 
the  perturbing  acceleration  due  to  the  ocean  tide  For  the  particular  example  chosen,  the 

two  equations  are  accurate  to  better  than  two  percent,  which  is  suspected  to  be  a 

coincidence.  But  we  are  certain  that  this  example  reveals  the  order  of  magnitude  of  a 

relationship  existing  between  the  gravitational  effects  of  the  ocean  tide  and  the  solid  earth 

tide. 
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